library(drc)
library(readxl)
library(SuperExactTest)
library(VennDiagram)
library(DT)
library(Glimma)
library(vioplot)
library(plotly)
load("~/Documents/RNAseq/RNAseq_S2/RNaseq_environnement.RData")
'%!in%' <- function(x,y)!('%in%'(x,y))

On va comparer les données d’expression des cellules S2 en présence de Shavenbaby Activateur et Shavenbaby Répresseur et les données dans les embryons mutant svb . Pour cela, je cherche:

1 Les données des gènes cibles dans les embryons

## Les données de toutes les amorces
data_embryon_all_probes = read_excel("~/Documents/Svb_Embryon/Data_Svbmutant_vs_WT_embryos_allProbes_Ced.xlsx")

#datatable(data_embryon_all_probes,rownames = F,filter = "top")

1.1 Les données des gènes différentiellement exprimés

data_embryon_selected = read_excel("~/Documents/Svb_Embryon/Data_Svbmutant_vs_WT_embryos_selected.xlsx")
## New names:
## * SYMBOL -> SYMBOL...3
## * SYMBOL -> SYMBOL...6
datatable(data_embryon_selected,rownames = F,filter = "top")

1.2 Les données des gènes cibles embryonnaires épidermiques (Menoret et al. 2013)

data_embryon_epiderm = read_excel("~/Documents/Svb_Embryon/Epidermal_embryonic_Svb_targets_updated.xlsx")
data_embryon_epiderm = as.data.frame(data_embryon_epiderm)

datatable(data_embryon_epiderm,rownames = F,filter = "top",caption =  "Genes Cibles Embryonnaires")
embryon_epiderm_tab = Dm6_gene[names(Dm6_gene)%in%data_embryon_epiderm$SYMBOL,]
embryon_epiderm_tab = as.data.frame(embryon_epiderm_tab)

#Obtention des groupes controles 
noembryon_epiderm_tab = Dm6_gene[names(Dm6_gene)%!in%data_embryon_epiderm$SYMBOL,]
noembryon_epiderm_tab = as.data.frame(noembryon_epiderm_tab)
noembryontab_random = noembryon_epiderm_tab[sample(1:nrow(noembryon_epiderm_tab),149),]

1.2.1 Les gènes cibles directes

1.2.1.1 8-10 h

gene_cibles_embryon_8_10 = data_embryon_epiderm$SYMBOL[data_embryon_epiderm$`ChIP peaks in 5kb (8-10h)` =="yes"]



embryon_epiderm_tab_8_10 = Dm6_gene[names(Dm6_gene)%in%gene_cibles_embryon_8_10,]
embryon_epiderm_tab_8_10 = as.data.frame(embryon_epiderm_tab_8_10)


datatable(as.data.frame(gene_cibles_embryon_8_10),rownames = F,filter = "top", caption = "Genes Cibles Embryon 8-10h")
notarget_set_tab_random_112 = noembryon_epiderm_tab[sample(1:nrow(noembryon_epiderm_tab),112),]

1.2.1.2 12-14h

gene_cibles_embryon_12_14 = data_embryon_epiderm$SYMBOL[data_embryon_epiderm$`ChIP peaks in 5kb (12-14h)`=="yes"]

embryon_epiderm_tab_8_10 = Dm6_gene[names(Dm6_gene)%in%gene_cibles_embryon_8_10,]
embryon_epiderm_tab_8_10 = as.data.frame(embryon_epiderm_tab_8_10)


datatable(as.data.frame(gene_cibles_embryon_12_14),rownames = F,filter = "top", caption = "Genes Cibles Directes Embryon 12-14h")
notarget_set_tab_random_85 = noembryon_epiderm_tab[sample(1:nrow(noembryon_epiderm_tab),85),]

2 Les données des gènes cibles dans les cellules S2

2.1 Les gènes cibles

Dans les cellules S2, un gène cible est activé par l’Activateur, réprimé par le Répresseur et l’effet antagoniste de Shavenbaby est significatif

upAct = row.names(lrt.2.tables$ActVsContr)[lrt.2.tables$ActVsContr$logFC > 0 ]
downRep =  row.names(lrt.2.tables$RepVsContr)[lrt.2.tables$RepVsContr$logFC <0 ]
signActRep = row.names(lrt.2.tables$ActVsRep)[lrt.2.tables$ActVsRep$BH.bin == 1]
t = upAct[upAct %in% downRep]
gene_cible_S2 = t[t%in%signActRep]
datatable(as.data.frame(gene_cible_S2),rownames = F,filter = "top", colnames = "Genes Cibles Cellules S2")

2.2 Les gènes cibles directes

Pour déterminer les gènes cibles directes de Shavenbaby, deux méthodes ont été utilisées, Soit en réattribuant le pic au gène avec ChIPpeakanno puis en faisant un superTestExact.

2.2.1 Gènes cibles avec la méthode de la distance ( coordonnées centriques)

Les gènes cibles sont des gènes dont la distance des pics et des gènes cibles directes est comprise entre 0 et 5kb

genes_cibles_S2_distance = read.table("/Users/alexmanchenoferris/Documents/ListeGene/target_with_pics/target_gene_with_peaks.txt", header = T)
datatable(genes_cibles_S2_distance[,3:4],rownames = F,filter = "top")
genes_cibles_S2_distance_direct = unique(genes_cibles_S2_distance$gene[genes_cibles_S2_distance$distance<=5000])

2.2.2 Gènes cibles directes (pics centriques)

genes_cibles_S2_supertest = read.csv("~/Documents/ListeGene/tableau_testMultiple/summary.table_upActdownRepActRepSign.csv",header = T)
genes_cibles_S2_supertest_direct = strsplit(genes_cibles_S2_supertest$Elements[genes_cibles_S2_supertest$Degree==4] ,", ")
datatable(as.data.frame(genes_cibles_S2_supertest_direct[[1]]),rownames = F,filter = "top", caption = "Genes cibles directes avec la méthode SuperExactTest")

3 Comparaison

3.1 Les gènes cibles S2

3.1.1 vs gènes différentiellement exprimés

venn.diagram(x = list(S2 = gene_cible_S2, embryo = data_embryon_selected$SYMBOL...3), main = "Comparaison Genes Cibles S2 et Embryon ", filename = "comparaison_cibles_S2embryo_tot.png",imagetype = "png", height = 5000, width = 5000, print.mode = c("raw","percent"), col = c("orange","blue"),fill = c("orange","blue"))
## [1] 1
commun_DE_select = gene_cible_S2[gene_cible_S2%in%data_embryon_selected$SYMBOL...3]
Comparaison des cibles S2 et embryon

Il y a 17 gènes communs entre les gènes cibles S2 et les gènes différentiellement exprimés en embryon.: Acer, CG1544, CG15818, CG17119, CG31436, CG32521, CG32700, CG7384, CG9119, Cyp6a20, Cyp6a8, GC1, Gip, GstD2, Karl, Orcokinin, Zip99C

3.1.2 vs gènes embryonnaires épidermiques (Menoret et al. 2013)

venn.diagram(x = list(S2 = gene_cible_S2, embryo = data_embryon_epiderm$SYMBOL), main = "Comparaison Genes Cibles  S2 et Embryon ", filename = "comparaison_cibles_S2embryo_Delphine.png",imagetype = "png", height = 5000, width = 5000, print.mode = c("raw","percent"), col = c("orange","blue"),fill = c("orange","blue"))
## [1] 1
commun_DE_epiderm = gene_cible_S2[gene_cible_S2%in%data_embryon_epiderm$SYMBOL]
Comparaison des cibles directes S2 et embryon

Il y a 3 gènes communs entre les gènes cibles S2 et les gènes épidermiques : CG15818, CG3831, CG8112

3.2 Les gènes cibles directes (coordonnées centriques)

3.2.1 vs gènes différentielles exprimés dans les embryons

venn.diagram(x = list(S2 = genes_cibles_S2_distance_direct, embryo = data_embryon_selected$SYMBOL...3), main = "Comparaison Genes Cibles Directes S2 (distance) et Embryon ", filename = "comparaison_ciblesdirecte_S2embryo_tot.png",imagetype = "png", height = 5000, width = 5000, print.mode = c("raw","percent"), col = c("orange","blue"),fill = c("orange","blue"))
## [1] 1
commun_DE_epiderm_distance = genes_cibles_S2_distance_direct[genes_cibles_S2_distance_direct%in%data_embryon_selected$SYMBOL...3]
Comparaison des cibles directes S2 (peaks centric) et gènes cibles embryonnaires

16 gènes différentiellement exprimés dans les embryons sont également gènes cibles directes de Shavenbaby: Acer, CG7384, Cyp6a20, Orcokinin, GC1, Zip99C, CG1544, Gip, Karl, CG32521, CG15818, Cyp6a8, GstD2, CG17119, CG31436, CG32700

3.2.2 vs gènes cibles épidermiques (Menoret et al. 2013)

venn.diagram(x = list(S2 = genes_cibles_S2_distance_direct, embryo = data_embryon_epiderm$SYMBOL), main = "Comparaison Genes Cibles Directes S2 (distance) et Embryon (épidermiques) ", filename = "comparaison_ciblesdirecte_S2embryo_epiderm.png",imagetype = "png", height = 5000, width = 5000, print.mode = c("raw","percent"), col = c("orange","blue"),fill = c("orange","blue"))
## [1] 1
commun_DE_epiderm_distance = genes_cibles_S2_distance_direct[genes_cibles_S2_distance_direct%in%data_embryon_epiderm$SYMBOL]
Comparaison des cibles directes S2 (peaks centric) et gènes cibles embryonnaires

3 gènes différentiellement exprimés dans les embryons sont également gènes cibles directes de Shavenbaby: CG3831, CG8112, CG15818

3.2.3 vs gènes cibles directes 8-10h

venn.diagram(x = list(S2 = genes_cibles_S2_distance_direct, embryo = gene_cibles_embryon_8_10), main = "Comparaison Genes Cibles Directes S2 et Embryon ( 8-10h) ", filename = "comparaison_ciblesdirecte_S2embryo_8_10.png",imagetype = "png", height = 5000, width = 5000, print.mode = c("raw","percent"), col = c("orange","blue"),fill = c("orange","blue"))
## [1] 1
Comparaison des cibles directes S2 et embryon (8-10h)

3.2.4 vs gènes cibles directes 12-14h

venn.diagram(x = list(S2 = genes_cibles_S2_distance_direct, embryo = gene_cibles_embryon_12_14), main = "Comparaison Genes Cibles Directes S2 et Embryon ( 12-14h) ", filename = "comparaison_ciblesdirecte_S2embryo_1214.png",imagetype = "png", height = 5000, width = 5000, print.mode = c("raw","percent"), col = c("orange","blue"),fill = c("orange","blue"))
## [1] 1
Comparaison des cibles directes S2 et embryon Delphine 12-14h
commun_8_10 = genes_cibles_S2_distance$gene[genes_cibles_S2_distance$gene%in%gene_cibles_embryon_8_10]



commun_12_14 = genes_cibles_S2_distance$gene[genes_cibles_S2_distance$gene%in%gene_cibles_embryon_12_14]

Les gènes cibles directes de Shavenbaby embryonnaires et des cellules S2 sont différents. Il y a 3 gènes communs CG3831, CG8112, CG15818 quand on compare les gènes cibles directs des S2 avec les gènes cibles directs embryonnaires 8-10h et 1 seul CG8112 avec les gènes cibles directs 12-14h.

3.2.5 Gènes cibles directes (pics centriques)

venn.diagram(x = list(S2 = genes_cibles_S2_supertest_direct[[1]], embryo = gene_cibles_embryon_8_10), main = "Comparaison Genes Cibles Directes S2 ( SuperTest) et Embryon ( 8-10h) ", filename = "comparaison_ciblesdirecte_S2embryo_supertest_8_10.png",imagetype = "png", height = 5000, width = 5000, print.mode = c("raw","percent"), col = c("orange","blue"),fill = c("orange","blue"))
## [1] 1
commun_8_10 = genes_cibles_S2_supertest_direct[[1]][genes_cibles_S2_supertest_direct[[1]]%in%gene_cibles_embryon_8_10]
venn.diagram(x = list(S2 = genes_cibles_S2_supertest_direct[[1]], embryo = gene_cibles_embryon_12_14), main = "Comparaison Genes Cibles Directes S2 (SuperTest) et Embryon ( 12-14h) ", filename = "comparaison_ciblesdirecte_S2embryosuperTest_1214.png",imagetype = "png", height = 5000, width = 5000, print.mode = c("raw","percent"), col = c("orange","blue"),fill = c("orange","blue"))
## [1] 1
commun_12_14 = genes_cibles_S2_supertest_direct[[1]][genes_cibles_S2_supertest_direct[[1]]%in%gene_cibles_embryon_12_14]

Comparaison des cibles directs S2 et embryon 8-10h

Comparaison des cibles directs S2 et embryon 12-14h

Les gènes cibles directes de Shavenbaby embryonnaires et des cellules S2 sont différents. Il y a 3 gènes communs CG15818, CG3831, CG8112 quand on compare les gènes cibles directs des S2 avec les gènes cibles directs embryonnaires 8-10h et 1 seul CG8112 avec les gènes cibles directs 12-14h.

4 Comparaison de la fixation de Shavenbaby sur les enhancers

Les enhancers embryonnaires sont ceux déterminés dans l’étude (Menoret et al. 2013) et ceux en S2 sont ceux de (Arnold et al. 2013) .

4.1 Sans nettoyage

4.1.1 Par bases

coverage_enhancer_embryon = read.table("/Users/alexmanchenoferris/Documents/coverage/coverage_S2_Embryon/Coverage_SvbS2_invivo_enhancerEmbryonnaire.tab")
colnames(coverage_enhancer_embryon) = c("chr","start","end","cov_SvbS2", "cov_SvbEmbryo")


fig = plot_ly(y = coverage_enhancer_embryon$cov_SvbS2, type = "box" , name = "Svb (S2)", boxpoints = FALSE)
fig = fig %>% add_trace(y = coverage_enhancer_embryon$cov_SvbEmbryo, name= "Svb (in vivo)", boxpoints = FALSE)
fig = fig %>% layout(title = "Coverage of Svb on embryonic enhancer")
fig
plot(x = coverage_enhancer_embryon$cov_SvbEmbryo, coverage_enhancer_embryon$cov_SvbS2)

coverage_enhancer_S2 = read.table("/Users/alexmanchenoferris/Documents/coverage/coverage_S2_Embryon/multiBamSummary/Coverage_SvbS2_invivo_enhancerS2.tab")
colnames(coverage_enhancer_S2) = c("chr","start","end","cov_SvbS2", "cov_SvbEmbryo")
fig = plot_ly(y = coverage_enhancer_S2$cov_SvbS2, type = "box" , name = "Svb (S2)", boxpoints = FALSE)
fig = fig %>% add_trace(y = coverage_enhancer_S2$cov_SvbEmbryo, name= "Svb (in vivo)", boxpoints = FALSE)
fig = fig %>% layout(title = "Coverage of Svb on S2 enhancer")
fig
plot(x = coverage_enhancer_S2$cov_SvbEmbryo, coverage_enhancer_S2$cov_SvbS2)

4.2 Par Régions

coverage_enhancer_embryon_region = read.table("/Users/alexmanchenoferris/Documents/coverage/coverage_S2_Embryon/multiBamSummary/Coverage_SvbS2_invivo_enhancerinvivo.tab", header = F)
colnames(coverage_enhancer_embryon_region) = c("chr","start","end","cov_SvbS2", "cov_SvbEmbryo")
enhancer_embryon = read.table("/Users/alexmanchenoferris/Documents/Svb_Embryon/bedfiles/enhancer_epidermique_dm6.bed")
enhancer_embryon = enhancer_embryon[order(enhancer_embryon$V1),]
coverage_enhancer_embryon_region$name = enhancer_embryon$V4
fig = plot_ly(coverage_enhancer_embryon_region, x = ~name , type = "bar",y = ~cov_SvbS2, name = "Svb S2")
fig = fig %>% add_trace(y = ~cov_SvbEmbryo, name = "Svb Embryon")
fig <- fig %>% layout(yaxis = list(title = 'Count'), barmode = 'group')
fig
coverage_enhancer_S2_region = read.table("/Users/alexmanchenoferris/Documents/coverage/coverage_S2_Embryon/multiBamSummary/Coverage_SvbS2_invivo_enhancerS2.tab", header = F)
colnames(coverage_enhancer_S2_region) =c("chr","start","end","cov_SvbS2", "cov_SvbEmbryo")
coverage_enhancer_S2_region$name = rownames(coverage_enhancer_S2_region)
fig2 = plot_ly(coverage_enhancer_S2_region, x = ~name , type = "bar",y = ~cov_SvbS2, name = "Svb S2")
fig2 = fig2 %>% add_trace(y = ~cov_SvbEmbryo, name = "Svb Embryon")
fig2 <- fig2 %>% layout(yaxis = list(title = 'Count'), barmode = 'group')
fig2

4.3 Nettoyant le bruit

Le Coverage est calculé sur les fichiers dont le bruit a été traité en normalisant sur la mappabilité du signal et en divisant IP/Input tout cela en log2.

4.3.1 Par Régions

coverage_enhancer_embryon_region_nettoye = read.table("/Users/alexmanchenoferris/Documents/coverage/coverage_S2_Embryon/log2IPInput/Coverage_svb_S2_embryo_refEmbryo.tab", header = F)
colnames(coverage_enhancer_embryon_region_nettoye) = c("chr","start","end","cov_SvbS2", "cov_SvbEmbryo")
enhancer_embryon = read.table("/Users/alexmanchenoferris/Documents/Svb_Embryon/bedfiles/enhancer_epidermique_dm6.bed")
enhancer_embryon = enhancer_embryon[order(enhancer_embryon$V1),]
coverage_enhancer_embryon_region_nettoye$name = enhancer_embryon$V4
fig = plot_ly(coverage_enhancer_embryon_region_nettoye, x = ~name , type = "bar",y = ~cov_SvbS2, name = "S2",color = I ( "green"))
fig = fig %>% add_trace(y = ~cov_SvbEmbryo, name = "Embryon", color = I("purple"))
fig <- fig %>% layout(yaxis = list(title = 'Svb Enrichment'), barmode = 'group', title = " Enrichment on Embryonic Enhancer")
fig

Pour les enhancers dans les S2, je n’en utilise que 42 pour avoir des groupes comparables entre les deux types cellulaires.

coverage_enhancer_S2_nettoye = read.table("/Users/alexmanchenoferris/Documents/coverage/coverage_S2_Embryon/log2IPInput/Coverage_svb_S2_embryo_refS2.tab")
colnames(coverage_enhancer_S2_nettoye) = c("chr","start","end","cov_SvbS2", "cov_SvbEmbryo")
coverage_enhancer_S2_nettoye$name = rownames(coverage_enhancer_S2_nettoye)

fig2 = plot_ly(coverage_enhancer_S2_nettoye, x = ~name , type = "bar",y = ~cov_SvbS2, name = "S2", color = I("green"))
fig2 = fig2 %>% add_trace(y = ~cov_SvbEmbryo, name = "Embryon",color = I("purple"))
fig2 <- fig2 %>% layout(yaxis = list(title = 'Svb Enrichment'), title = "Enrichment on S2 enhancers ",barmode = 'group')
fig2

Nous observons que sur les enhancers embryonnaires le signal en Svb est plus fort en embryon qu’en S2 Ce qui confirme notre observable que Svb se fixe moins en S2 sur les enhancers embryonnaires qu’en embryon.

Cette observation est inversée dans le cas des enhancers en S2 # References {.tabset .tabset-fade .tabset-pills}

Arnold, Cosmas D, Daniel Gerlach, Christoph Stelzer, Łukasz M Boryń, Martina Rath, and Alexander Stark. 2013. “Genome-wide quantitative enhancer activity maps identified by STARR-seq.” Science (New York, N.Y.) 339 (6123): 1074–7. https://doi.org/10.1126/science.1232542.

Menoret, Delphine, Marc Santolini, Isabelle Fernandes, Rebecca Spokony, Jennifer Zanet, Ignacio Gonzalez, Yvan Latapie, et al. 2013. “Genome-wide analyses of Shavenbaby target genes reveals distinct features of enhancer organization.” Genome Biology 14 (8): R86. https://doi.org/10.1186/gb-2013-14-8-r86.